#include "HY1C_out.h"
/* ========================================================================
 *
 *   MSl1b2info input-l1b-filename
 *
 * Description:
 * 
 * Modification history:
 *
 *     Programmer     Organization      Date       Description of change
 *   --------------   ------------    --------     ---------------------
 *   Bryan A. Franz   GSC             28 March 1999 Original development
 *   Dan Knowles      SAIC            2004          Complete rewrite to create geoboxes
 *   W. Robinson      SAIC            10 dec 2004   add CZCS sensor
 *   Dan Knowles      SAIC            18 Aug 2005   Changed daynight logic to be based on solar zenith angle
 *
 * ======================================================================== */

#include <stdio.h>
#include <sys/types.h>
#include <sys/stat.h>
#include <stdlib.h>
#include <fcntl.h>
#include <string.h>
#include <errno.h>
#include <unistd.h>
#include <libgen.h>
#include <math.h>

#include "l12_proto.h"
#include "time_utils.h"

#define NUM_SPACECRAFT_DIRECTION_MODES 3
#define DEFAULT_SPACECRAFT_DIRECTION 0   /* set to an invalid value */
#define ASCENDING_TRACK  1
#define DESCENDING_TRACK 2

#define NORTH 1
#define SOUTH 0

#ifdef TRUE
#undef TRUE
#endif

#define TRUE 1
#define FALSE 0

#define TMP_FILENAME_MAX 255

#define DEFAULT_COORD_VALUE -999.   /* set to an invalid value */

#define DEFAULT_DAYNIGHT_MODE 0
#define DAY_MODE (1<<0)
#define NIGHT_MODE (1<<1)

/* equatorial radius, didn't bother to use polar value also since precision not that important for current use */
#define EARTH_RADIUS_EQUATORIAL 6378   


typedef struct {
    float32 north_lat;
    float32 south_lat;
    float32 west_lon;
    float32 east_lon;
    unsigned char daynightflag;
    uint32_t pixel_count;
} box_t;

typedef struct {
    int sensor_id;
    int day_node;
} day_node_t;


void set_north_south_boundaries(float32, float32 *, float32 *);
void set_west_east_boundaries(float32, float32 *, float32 *);
int check_if_in_west_east_boundaries(float32, float32, float32);
void get_box_square_km(double, double, double, double, double *);
void get_box_square_km_sublevel(double, double, double, double, double *);
int check_if_in_box(float32, float32, float32, float32, float32, float32);
double get_lon_distance(double,double);

#define USAGESTR \
"%s %s (%s %s)\n\n\
Usage: %s [-n number-of-boxes] [-s] L1-filename [met-file]\n\n\
Where:\n\
\t-d degree limits to set on boxes (use instead of -n)\n\
\t-i set subsampling increment (default=1 gets all scan lines and pixels)\n\
\t-n generates up to this number of lat/lon boxes per orbit-side\n\
\t-s prints only data needed for database fields\n\
\t-v verify number of pixels that match in boxed region\n\
\n\
"


#define PRINT_USAGE(x)  HY1C_out(USAGESTR,(x),VERSION,__DATE__,__TIME__,(x))


  

/* -------------------------------------------------------------------- *
 *                              main                                    *
 * -------------------------------------------------------------------- */
int main(int argc, char *argv[])
{
    l1str       l1rec;              /* generic level-1b scan structure      */
    filehandle  l1file;             /* input file handle                    */
    instr       input;
    float64     utime;
    char       *ztime;
#if 0
    char        geofile[FILENAME_MAX] = "";   
    char        geofile[TMP_FILENAME_MAX] = "";
#endif
    int16       mm, dd;

    int32_t   spix     = 0;          /* start pixel for subscene process   */
    int32_t   cpix     = -1;         /* center pixel of subscene scan      */
    int32_t   epix     = -1;         /* end pixel for subscene process     */
    int32_t   sscan    = 0;          /* start scan for subscene process    */
    int32_t   cscan    = -1;         /* center scan for subscene process   */
    int32_t   escan    = -1;         /* end scan for subscene process      */
    int32_t   iscan;
    int32_t   ip;
    int32_t  num_match = 0;
    int32_t  num_no_match = 0;
    int32_t  num_nav_fail = 0;
    int32_t  num_bad_lat_lon_pixels = 0;
    int32_t  num_pixels = 0;

    int   i;
    int   box_index;
    int   spacecraft_direction_index;
    int   curr_spacecraft_direction = DEFAULT_SPACECRAFT_DIRECTION;
    int   prev_spacecraft_direction = DEFAULT_SPACECRAFT_DIRECTION;
    int   initial_spacecraft_direction = DEFAULT_SPACECRAFT_DIRECTION;
    int   final_spacecraft_direction = DEFAULT_SPACECRAFT_DIRECTION;
    int   match;
    int   in_box;
    int   num_boxes = 0;
    int   num_degrees = 0;
    int   show_diagnostic_info;
    int   show_standard_info;
    int   show_sdps_info;
    int   case_n = FALSE;
    int   case_s = FALSE;
    int   case_d = FALSE;
    int   case_v = FALSE;
    int   total_box_count = 0;
    int   user_defined_increment = 1;
    int   exitflag = 0;

    double total_box_square_km = 0;
    double box_square_km;
    double percent_match;

    unsigned char curr_daynight_mode = DEFAULT_DAYNIGHT_MODE;
    unsigned char master_daynightflag = DEFAULT_DAYNIGHT_MODE;

    static char *node[3] = {"Unidentified","Ascending","Descending"};
    static char *daynight_string[4] = {"Unidentified","Day","Night","Both"};

    float32 epix_northern_lat = DEFAULT_COORD_VALUE;
    float32 epix_southern_lat = DEFAULT_COORD_VALUE;
    float32 spix_northern_lat = DEFAULT_COORD_VALUE;
    float32 spix_southern_lat = DEFAULT_COORD_VALUE;
    float32 northern_lat = DEFAULT_COORD_VALUE;
    float32 southern_lat = DEFAULT_COORD_VALUE;
    float32 eastern_lon  = DEFAULT_COORD_VALUE;
    float32 western_lon  = DEFAULT_COORD_VALUE;
    float32 center_lon   = DEFAULT_COORD_VALUE;
    float32 center_lat   = DEFAULT_COORD_VALUE;
    float32 prev_lat_cpix = DEFAULT_COORD_VALUE;
    float32 tmp_north_lat;
    float32 tmp_south_lat;
    float32 tmp_west_lon;
    float32 tmp_east_lon;
    float32 lat_breakpoint;
    float32 increment;


    box_t **box = NULL;  /* will become 2 dimensional box array in [boxes][direction]  */
    int counter;
    int num_args;

    char *option_string = "d:hi:n:sv";
    int options =0;
    while ((options = getopt(argc, argv, option_string)) != -1)
    {
        switch (options)
        {
            case 'd':
                case_d = TRUE;
                num_degrees = atoi(optarg);
                break;
            case 'i':
                user_defined_increment = atoi(optarg);
                break;
            case 'h':
                PRINT_USAGE(argv[0]);
                exit(FATAL_ERROR);
                break;
            case 'n':
                case_n = TRUE;
                num_boxes = atoi(optarg);
                break;
            case 's':
                case_s = TRUE;
                break;
            case 'v':
                case_v = TRUE;
                break;
            default:
                break;
        }
    }                          


    /*******************************************************************************************
    *    Check options and set num_boxes
    *******************************************************************************************/

    if (case_d == TRUE)
    {
        if (num_degrees >= 1 && num_degrees <= 180)
        {
            if (case_n != TRUE)
            {
                num_boxes = (int) 180 / num_degrees;
                HY1C_out("num_boxes=%d\n",num_boxes);
            }
            else
            {
                HY1C_out("INPUT ERROR: Cannot choose both -d and -n options\n");
                PRINT_USAGE(argv[0]);
                exit(FATAL_ERROR);
            }
        }
        else
        {
            HY1C_out("INPUT ERROR: -d (degrees lat) option must be between 1 and 180\n");
            PRINT_USAGE(argv[0]);
            exit(FATAL_ERROR);
        }
    } 
    else if (case_n == TRUE)
    {
        if (num_boxes < 0)
        {
            HY1C_out("INPUT ERROR: -n (number of lat divisions) option must be greater than or equal to 0\n");
            PRINT_USAGE(argv[0]);
            exit(FATAL_ERROR);
        }
    }
    else
    {
        num_boxes = 0;
    }
/*    else if (case_s == TRUE)
    {
        if (case_n != TRUE && case_d != TRUE)
        {
            HY1C_out("INPUT ERROR: -s option requires -d or -n option also\n");
            PRINT_USAGE(argv[0]);
            exit(FATAL_ERROR);
        }
    }
*/    

    /*******************************************************************************************
    *    Setup print options
    *******************************************************************************************/

    if (case_s == TRUE)
    {
        show_sdps_info = TRUE;
        show_standard_info = FALSE;
        show_diagnostic_info = FALSE;
    } 
    else if (case_n == TRUE || case_d == TRUE)
    {
        show_sdps_info = FALSE;
        show_standard_info = TRUE;
        show_diagnostic_info = TRUE;
    } 
    else
    {
        show_sdps_info = FALSE;
        show_standard_info = TRUE;
        show_diagnostic_info = FALSE;
    }
        

    num_args = argc - optind; 

    msl12_input_init(&input);
    filehandle_init(&l1file);
    l1file.input = &input;

    /*									*/
    /* Process required command-line arguments                          */
    /*									*/
    switch (num_args) {
      case 2:
        strcpy(l1file.name,argv[optind+0]);
        /* strcpy(geofile,argv[optind+1]);*/
        l1file.geofile = argv[optind+1];
        break;
      case 1:
        strcpy(l1file.name,argv[optind+0]);
        break;
      default:
        PRINT_USAGE(argv[0]);
        exit(FATAL_ERROR);
        break;
    }

    /* Set default input parameters */
    if (msl12_input_defaults(&l1file, &input) != 0) {
        HY1C_out("-E- %s: Error parsing input parameters.\n",argv[0]);
        exit(FATAL_ERROR);
    }

    input.proc_sst=0;
    input.proc_ocean=0;
    input.atmocor=0;
    input.proc_land=0;

    /*									*/
    /* Open input file and get sensor and scan information from handle  */
    /*                                                                  */
    /* done in msl12_input_defaults()
    if (getFormat(&l1file) == -1) {
        HY1C_out("-E- %s: Unrecognized input file format.\n",argv[0]);
        exit(FATAL_ERROR);
    }
    */
    if (openl1(&l1file) != 0) {
        HY1C_out("-E- %s: Error opening %s for reading.\n",argv[0],l1file.name);
        exit(FATAL_ERROR);
    }

    /*									*/
    /* Allocate memory for L1 scan data			  	        */
    /*									*/
    if ( alloc_l1(l1file.npix,NBANDS,NBANDSIR,l1file.n_inprods,&l1rec) == 0 ) {
        HY1C_out("-E- %s: Unable to allocate L1 record.\n",argv[0]);
        exit( FATAL_ERROR );
    }

    spix = 0;
    epix = l1file.npix-1;
    cpix = spix  + (epix-spix+1)/2;

    sscan = 0;
    escan = l1file.nscan-1;
    cscan = sscan + (escan-sscan+1)/2;

    l1file.spix = spix;
    l1file.epix = epix;
    l1file.dpix = 1;

    if (show_standard_info == TRUE)
    {
        HY1C_out("Sensor=%s\n",sensorName[l1file.sensorID]);
    }

    if (show_sdps_info == TRUE || show_standard_info == TRUE)
    {
        HY1C_out("Number_of_Scans=%d\n",l1file.nscan);
    }

    if (show_standard_info == TRUE)
    {
        HY1C_out("Number_of_Pixels_per_Scan=%d\n",l1file.npix);
        /*
        HY1C_out("Percent_Clouds=%d\n",l1file.percent_cloud);
        HY1C_out("Percent_Land=%d\n",l1file.percent_land);
        HY1C_out("Percent_Water=%d\n",l1file.percent_water);
        */
    }


    if (num_boxes > 0)
    {
        box = (box_t **) calloc(num_boxes,sizeof(box_t*));

        if(box == NULL)
        {
            HY1C_out("out of memory\n");
            exit(FATAL_ERROR);
        }
    }

    for(counter=0;counter<num_boxes;counter++)
    {
        box[counter] = calloc(NUM_SPACECRAFT_DIRECTION_MODES,sizeof(box_t));
        if(box[counter] == NULL)
        {
            HY1C_out("out of memory\n");
            exit(FATAL_ERROR);
        }
    }



    /***********************************************************************************************
    *    Initialize boxes
    ***********************************************************************************************/

    for (box_index = 0; box_index < num_boxes; box_index++)
    {
        for (spacecraft_direction_index = 0; spacecraft_direction_index < NUM_SPACECRAFT_DIRECTION_MODES; spacecraft_direction_index++)
        {
            box[box_index][spacecraft_direction_index].north_lat = DEFAULT_COORD_VALUE;
            box[box_index][spacecraft_direction_index].south_lat = DEFAULT_COORD_VALUE;
            box[box_index][spacecraft_direction_index].west_lon = DEFAULT_COORD_VALUE;
            box[box_index][spacecraft_direction_index].east_lon = DEFAULT_COORD_VALUE;
            box[box_index][spacecraft_direction_index].daynightflag = DEFAULT_DAYNIGHT_MODE;
        }
    }


    /***************************************************************************
    *    Loop through all scan lines
    ***************************************************************************/

/*
    for (iscan=sscan; iscan<=escan; iscan++) 
    {  */

    iscan=sscan;

    while (iscan<=escan)
    {
        readl1(&l1file,iscan,&l1rec);

        if (iscan == sscan) {
            utime = yds2unix(*l1rec.year,*l1rec.day,*l1rec.msec/((double) 1000.0));
            ztime = ydhmsf(utime,'G');
            yd2md(*l1rec.year,*l1rec.day,&mm,&dd);
            if (show_standard_info == TRUE)
            {
                HY1C_out("\n");
                HY1C_out("INFO: FIRST SCAN LINE\n");
                HY1C_out("Start_Date=%4.4d-%2.2d-%2.2d %.2s:%.2s:%.2s.%.3s\n",
                        *l1rec.year,mm,dd,ztime+7,ztime+9,ztime+11,ztime+13);
                HY1C_out("Start_Year=%.4s\n",ztime);
                HY1C_out("Start_Day=%.3s\n",ztime+4);
                HY1C_out("Start_Hour=%.2s\n",ztime+7);
                HY1C_out("Start_Min=%.2s\n",ztime+9);
                HY1C_out("Start_Sec=%.2s\n",ztime+11);
                HY1C_out("Start_Msec=%.3s\n",ztime+13);
                HY1C_out("Start_Time=%s\n",ztime);
                HY1C_out("Upper_Left_Lon=%f\n",l1rec.lon[spix]);
                HY1C_out("Upper_Left_Lat=%f\n",l1rec.lat[spix]);
                HY1C_out("Upper_Right_Lon=%f\n",l1rec.lon[epix]);
                HY1C_out("Upper_Right_Lat=%f\n",l1rec.lat[epix]);
            }

            if (show_sdps_info == TRUE)
            {
                HY1C_out("Start_Date=%4.4d-%2.2d-%2.2d\n",*l1rec.year,mm,dd);
                HY1C_out("Start_Time=%.2s:%.2s:%.2s.%.3s\n",ztime+7,ztime+9,ztime+11,ztime+13);
            }

        }
    
        if (iscan == escan) {
            utime = yds2unix(*l1rec.year,*l1rec.day,*l1rec.msec/((double) 1000.0));
            ztime = ydhmsf(utime,'G');
            yd2md(*l1rec.year,*l1rec.day,&mm,&dd);

            if (show_standard_info == TRUE)
            {
                HY1C_out("\n");
                HY1C_out("INFO: LAST SCAN LINE\n");
                HY1C_out("End_Date=%4.4d-%2.2d-%2.2d %.2s:%.2s:%.2s.%.3s\n",
		       *l1rec.year,mm,dd,ztime+7,ztime+9,ztime+11,ztime+13);
                HY1C_out("End_Year=%.4s\n",ztime);
                HY1C_out("End_Day=%.3s\n",ztime+4);
                HY1C_out("End_Hour=%.2s\n",ztime+7);
                HY1C_out("End_Min=%.2s\n",ztime+9);
                HY1C_out("End_Sec=%.2s\n",ztime+11);
                HY1C_out("End_Msec=%.3s\n",ztime+13);
                HY1C_out("End_Time=%s\n",ztime);
                HY1C_out("Lower_Left_Lon=%f\n",l1rec.lon[spix]);
                HY1C_out("Lower_Left_Lat=%f\n",l1rec.lat[spix]);
                HY1C_out("Lower_Right_Lon=%f\n",l1rec.lon[epix]);
                HY1C_out("Lower_Right_Lat=%f\n",l1rec.lat[epix]);
	    }

            if (show_sdps_info == TRUE)
            {
                HY1C_out("End_Date=%4.4d-%2.2d-%2.2d\n",*l1rec.year,mm,dd);
                HY1C_out("End_Time=%.2s:%.2s:%.2s.%.3s\n",ztime+7,ztime+9,ztime+11,ztime+13);
            }

        }


        /***************************************************************************************************
        *    Determine MAX and MIN values for epix_lat and spix_lat based on current and all preceeding scan lines
        ***************************************************************************************************/

        
        if ((l1rec.flags[spix] & NAVFAIL) == 0)
        {
            set_north_south_boundaries(l1rec.lat[spix],&spix_northern_lat,&spix_southern_lat);
        }

        if ((l1rec.flags[epix] & NAVFAIL) == 0)
        {
            set_north_south_boundaries(l1rec.lat[epix],&epix_northern_lat,&epix_southern_lat);
        }

        /**************************************************************
        *    Find lat and lon of the center pixel of the middle scan
        **************************************************************/

        if (center_lat == DEFAULT_COORD_VALUE || center_lon == DEFAULT_COORD_VALUE)
        {
            if (iscan >= cscan && (l1rec.flags[cpix] & NAVFAIL) == 0 )
            {
                center_lat = l1rec.lat[cpix];
                center_lon = l1rec.lon[cpix];
            }
        }


        /***************************************************************************
        *    Determine spacecraft direction based on center pixel of the scan
        ***************************************************************************/

        if ((l1rec.flags[cpix] & NAVFAIL) == 0)
        {
            if (prev_lat_cpix != DEFAULT_COORD_VALUE)
            {
                if (l1rec.lat[cpix] > prev_lat_cpix)
                {
                    curr_spacecraft_direction = ASCENDING_TRACK;
                }
                else if (l1rec.lat[cpix] < prev_lat_cpix)
                {

                    curr_spacecraft_direction = DESCENDING_TRACK;
                }
                else
                {
                    curr_spacecraft_direction = prev_spacecraft_direction;
                }

                if (initial_spacecraft_direction == DEFAULT_SPACECRAFT_DIRECTION &&
                    curr_spacecraft_direction != DEFAULT_SPACECRAFT_DIRECTION)
                {
                    initial_spacecraft_direction = curr_spacecraft_direction;
                }

            }
            prev_lat_cpix  = l1rec.lat[cpix];
            prev_spacecraft_direction = curr_spacecraft_direction;

        }


        /***************************************************************************
        *    Loop through all pixels 
        ***************************************************************************/

        /* for (ip=spix; ip<=epix; ip++) */

        ip = spix;

        while (ip <= epix)
        {

            if ((l1rec.flags[ip] & NAVFAIL) == 0)
            {
                /***************************************************************************
                *    Determine day/night mode of pixel 
                ***************************************************************************/

                if (l1rec.solz[ip] < 90)
                {
                    curr_daynight_mode = DAY_MODE;
                }
                else
                {
                    curr_daynight_mode = NIGHT_MODE;
                }

                master_daynightflag |= curr_daynight_mode;


                /***************************************************************************************************
                *    Determine MAX and MIN values for both lat and lon based on current and all preceeding scan lines 
                ***************************************************************************************************/
        
                set_north_south_boundaries(l1rec.lat[ip],&northern_lat,&southern_lat);
                set_west_east_boundaries(l1rec.lon[ip],&western_lon,&eastern_lon);


                if (num_boxes > 0)
                {
                    /***************************************************************************
                    *    Determine and set boundaries for regular boxes
                    ***************************************************************************/
    
                    if (l1rec.lat[ip] >= -90. && l1rec.lat[ip] <= 90. && l1rec.lon[ip] >= -180. && l1rec.lon[ip] <= 180.)
                    {
                        i = 0; 
                        box_index = num_boxes; 
                        increment = 180. / num_boxes;
    
                        for (lat_breakpoint = 90. - increment; lat_breakpoint >= -90.; lat_breakpoint -= increment)
                        { 
                            if (l1rec.lat[ip] >= lat_breakpoint && l1rec.lat[ip] < (lat_breakpoint + increment) )
                            {
                                box_index = i;
                            }
                            i++;
    
                        }
 
                        if (box_index < num_boxes)
                        {
                            set_north_south_boundaries(l1rec.lat[ip],
                                                       &box[box_index][curr_spacecraft_direction].north_lat,
                                                       &box[box_index][curr_spacecraft_direction].south_lat);
                            set_west_east_boundaries(l1rec.lon[ip],
                                                     &box[box_index][curr_spacecraft_direction].west_lon,
                                                     &box[box_index][curr_spacecraft_direction].east_lon);
                            box[box_index][curr_spacecraft_direction].daynightflag |= curr_daynight_mode;
                            box[box_index][curr_spacecraft_direction].pixel_count++;
                        }
                    }
                    else
                    {
                        num_bad_lat_lon_pixels++;

                        if (l1rec.lat[ip] < -90. || l1rec.lat[ip] > 90.)
                        {
			    fHY1C_out(stderr, "WARNING_MSG: lat pixel out of range: lat=%f, pixel_index=%d", l1rec.lat[ip], ip);
                        }

                        if (l1rec.lon[ip] < -180. || l1rec.lon[ip] > 180.)
                        {
			    fHY1C_out(stderr, "WARNING_MSG: lon pixel out of range: lon=%f, pixel_index=%d", l1rec.lon[ip], ip);
                        }
                    }
                }
            }
            else
            {
                num_nav_fail++;    /* count failed pixels  */
            }

            num_pixels++;   /* count total pixels */

            if (ip < epix - user_defined_increment)
            {
                ip += user_defined_increment;  
            }
            else
            {
                ip++;
            }
        }


        if (iscan < escan - user_defined_increment)
        {
            iscan += user_defined_increment; 
        }
        else
        {
            iscan++;
        }
    }

    final_spacecraft_direction = curr_spacecraft_direction;


    if (num_bad_lat_lon_pixels > 0)
    {
        fHY1C_out(stderr, "WARNING_MSG: found %d lat/lon pixel(s) which were out of range", num_bad_lat_lon_pixels);
    }



    /***************************************************************************************************
    *    Merge initial scans box with appropriate box now that we know the initial_spacecraft_direction 
    ***************************************************************************************************/

    for (box_index = 0; box_index < num_boxes; box_index++)
    {
        if (box[box_index][DEFAULT_SPACECRAFT_DIRECTION].north_lat != DEFAULT_COORD_VALUE)
        {
            set_north_south_boundaries(box[box_index][DEFAULT_SPACECRAFT_DIRECTION].north_lat,
                                       &box[box_index][initial_spacecraft_direction].north_lat,
                                       &box[box_index][initial_spacecraft_direction].south_lat);
            set_north_south_boundaries(box[box_index][DEFAULT_SPACECRAFT_DIRECTION].south_lat,
                                       &box[box_index][initial_spacecraft_direction].north_lat,
                                       &box[box_index][initial_spacecraft_direction].south_lat);
            set_west_east_boundaries(box[box_index][DEFAULT_SPACECRAFT_DIRECTION].west_lon,
                                       &box[box_index][initial_spacecraft_direction].west_lon,
                                       &box[box_index][initial_spacecraft_direction].east_lon);
            set_west_east_boundaries(box[box_index][DEFAULT_SPACECRAFT_DIRECTION].east_lon,
                                       &box[box_index][initial_spacecraft_direction].west_lon,
                                       &box[box_index][initial_spacecraft_direction].east_lon);
            box[box_index][initial_spacecraft_direction].pixel_count += box[box_index][DEFAULT_SPACECRAFT_DIRECTION].pixel_count;
            box[box_index][initial_spacecraft_direction].daynightflag |= box[box_index][DEFAULT_SPACECRAFT_DIRECTION].daynightflag;

            /* re initialize default spacecraft direction box */
            box[box_index][DEFAULT_SPACECRAFT_DIRECTION].north_lat = DEFAULT_COORD_VALUE;
            box[box_index][DEFAULT_SPACECRAFT_DIRECTION].south_lat = DEFAULT_COORD_VALUE;
            box[box_index][DEFAULT_SPACECRAFT_DIRECTION].west_lon = DEFAULT_COORD_VALUE;
            box[box_index][DEFAULT_SPACECRAFT_DIRECTION].east_lon = DEFAULT_COORD_VALUE;
            box[box_index][DEFAULT_SPACECRAFT_DIRECTION].daynightflag = DEFAULT_DAYNIGHT_MODE;
        }
    }

    if (num_boxes != 0)
    {
        if (num_pixels == num_nav_fail)
        {
            fHY1C_out(stderr, "ERROR_MSG: all %d pixels have navigation failure\n", num_pixels);
            exitflag = 100;
        }
    
        else if (northern_lat == DEFAULT_COORD_VALUE ||
                 southern_lat == DEFAULT_COORD_VALUE ||
                 eastern_lon == DEFAULT_COORD_VALUE ||
                 western_lon == DEFAULT_COORD_VALUE)
        {
            fHY1C_out(stderr, "ERROR_MSG: lat/lon coordinates never set in this program\n");
            exitflag = 101;
        }

        else if (case_n == TRUE && (initial_spacecraft_direction == DEFAULT_SPACECRAFT_DIRECTION ||
                               final_spacecraft_direction == DEFAULT_SPACECRAFT_DIRECTION)
                )
        {
            fHY1C_out(stderr, "ERROR_MSG: initial/final direction never set in this program\n");
            exitflag = 102;
        }
    
    }

    /***************************************************************************************************
    *    Print summary info 
    ***************************************************************************************************/

    if (show_standard_info == TRUE)
    {
        HY1C_out("\n");
        HY1C_out("SUMMARY STATS\n");
        HY1C_out("Northernmost_Lat=%f\n",northern_lat);
        HY1C_out("Southernmost_Lat=%f\n",southern_lat);
        HY1C_out("Easternmost_Lon=%f\n",eastern_lon);
        HY1C_out("Westernmost_Lon=%f\n",western_lon);
        HY1C_out("Center_Lat=%f\n",center_lat);
        HY1C_out("Center_Lon=%f\n",center_lon);
        HY1C_out("Start_Node=%s\n",node[initial_spacecraft_direction]);
        HY1C_out("End_Node=%s\n",node[final_spacecraft_direction]);
        HY1C_out("Daynight=%s\n",daynight_string[master_daynightflag]);
        HY1C_out("Moon_in_SV=%1d\n",l1file.sv_with_moon);
    }


    if (show_diagnostic_info == TRUE)
    {
        HY1C_out("Total_Num_Pixels=%d\n",num_pixels);
        HY1C_out("Num_Nav_Fail_Pixels=%d\n",num_nav_fail);
        HY1C_out("\n");
        HY1C_out("End_of_Scan_Northernmost_Lat=%f\n",epix_northern_lat);
        HY1C_out("End_of_Scan_Southernmost_Lat=%f\n",epix_southern_lat);
        HY1C_out("Start_of_Scan_Northernmost_Lat=%f\n",spix_northern_lat);
        HY1C_out("Start_of_Scan_Southernmost_Lat=%f\n",spix_southern_lat);
        HY1C_out("\n");
    }
    


    if (show_sdps_info == TRUE)
    {
        HY1C_out("DayNightFlag=%d\n",master_daynightflag);
        HY1C_out("Moon_in_SV=%1d\n",l1file.sv_with_moon);
    }

    if (exitflag == 0 && master_daynightflag == 0)
    {
        exitflag = 103;
        fHY1C_out(stderr, "ERROR_MSG: master_daynightflag never set in this program\n");
    }


    /***************************************************************************************************
    *    Print box info 
    ***************************************************************************************************/

    for (spacecraft_direction_index = 0; spacecraft_direction_index < NUM_SPACECRAFT_DIRECTION_MODES; spacecraft_direction_index++)
    {
        for (box_index = 0; box_index < num_boxes; box_index++)
        {
            if (box[box_index][spacecraft_direction_index].north_lat != DEFAULT_COORD_VALUE)
            {
                total_box_count++;

                tmp_north_lat = box[box_index][spacecraft_direction_index].north_lat;
                tmp_south_lat = box[box_index][spacecraft_direction_index].south_lat;
                tmp_west_lon = box[box_index][spacecraft_direction_index].west_lon;
                tmp_east_lon = box[box_index][spacecraft_direction_index].east_lon;
              
                get_box_square_km((double) tmp_north_lat,
                                     (double) tmp_south_lat,
                                     (double) tmp_west_lon,
                                     (double) tmp_east_lon,
                                     &box_square_km);
 
                total_box_square_km += box_square_km;

                if (show_diagnostic_info == TRUE)
                {
                    HY1C_out("Box_Summary\n");
                    HY1C_out("Box_Northern_boundary=%f\n",box[box_index][spacecraft_direction_index].north_lat);
                    HY1C_out("Box_Southern_boundary=%f\n",box[box_index][spacecraft_direction_index].south_lat);
                    HY1C_out("Box_Western_boundary=%f\n",box[box_index][spacecraft_direction_index].west_lon);
                    HY1C_out("Box_Eastern_boundary=%f\n",box[box_index][spacecraft_direction_index].east_lon);
                    HY1C_out("Box_Daynightflag=%s\n",daynight_string[box[box_index][spacecraft_direction_index].daynightflag]);
                    HY1C_out("Box_Pixel_count=%d\n",box[box_index][spacecraft_direction_index].pixel_count);
                    HY1C_out("Box_size_in_square_million_km=%.2f\n",box_square_km/1000000.);
                    HY1C_out("\n");
                }
                if (show_sdps_info == TRUE)
                {
                    HY1C_out("GeoBox_%d=%f,%f,%f,%f,%d\n",
                           total_box_count,
                           box[box_index][spacecraft_direction_index].north_lat,
                           box[box_index][spacecraft_direction_index].south_lat,
                           box[box_index][spacecraft_direction_index].west_lon,
                           box[box_index][spacecraft_direction_index].east_lon,
                           box[box_index][spacecraft_direction_index].daynightflag);
                }
                
                
            }
        }
    }
    
    



    
    if (show_diagnostic_info == TRUE)
    {
        HY1C_out("Total_box_count=%d\n",total_box_count);
        HY1C_out("Total_box_size_square_million_km=%.2f\n",total_box_square_km/1000000.);
    }

    if (exitflag == 0 && case_n == TRUE && total_box_count == 0 && num_boxes != 0)
    {
        exitflag = 104;
        fHY1C_out(stderr, "ERROR_MSG: total_box_count = 0\n");
    }
  

    if (case_v == TRUE)
    {
 
        for (iscan=sscan; iscan<=escan; iscan++)
        {
            readl1(&l1file,iscan,&l1rec);
      
            for (ip=spix; ip<=epix; ip++)
            {
                if ( (l1rec.flags[ip] & NAVFAIL) == 0)
                {
                    match = FALSE;
        
                    for (box_index = 0; box_index < num_boxes; box_index++)
                    {
                        for (spacecraft_direction_index = 0; spacecraft_direction_index < NUM_SPACECRAFT_DIRECTION_MODES; spacecraft_direction_index++)
                        {
                            in_box = check_if_in_box(l1rec.lat[ip],
                                                     l1rec.lon[ip],
                                                     box[box_index][spacecraft_direction_index].north_lat,
                                                     box[box_index][spacecraft_direction_index].south_lat,
                                                     box[box_index][spacecraft_direction_index].west_lon,
                                                     box[box_index][spacecraft_direction_index].east_lon);
                            if (in_box == TRUE)
                            {
                                match = TRUE;
                            }
                        }
                    }
        
                    if (match == TRUE)
                    {
                        num_match++;
                    }
                    else
                    {
                        num_no_match++;
                    }
     
                }
            }
        }
  
        if ((num_match + num_no_match) > 0)
        {
            percent_match = 100. * num_match / (num_match + num_no_match); 
        }
        else
        {
            percent_match = 0; 
        }
 
        HY1C_out("percent_match=%f, num_match=%d, num_no_match=%d\n", percent_match, num_match, num_no_match);
    }


    /* free the dynamically allocated memory */

    if (num_boxes > 0)
    {
        for(counter = 0; counter < num_boxes; counter++)
        {
            free(box[counter]);
        }
        /* and finally free ptr to ptr....the space to store ptr to box_t*/
        free(box);
    }

    exit(exitflag);
}





/* take first arg (lat) and adjust northern_boundary or southern_boundary if needed  */

void
set_north_south_boundaries(float32 lat, float32 *northern_boundary, float32 *southern_boundary)
{
    if (lat > 90.)
    {
        lat = 90.;
    }
    
    if (lat < -90.)
    {
        lat = -90.;
    }

    if (*northern_boundary != DEFAULT_COORD_VALUE)
    {
        *northern_boundary = MAX(lat,*northern_boundary);
    }
    else
    {
        *northern_boundary = lat;
    }

    if (*southern_boundary != DEFAULT_COORD_VALUE)
    {
        *southern_boundary = MIN(lat,*southern_boundary);
    }
    else
    {
        *southern_boundary = lat;
    }

}

double 
get_lon_distance(double lon1,double lon2)
{
    double results;
    double distance_1;
    double distance_2;

    distance_1 = lon1 - lon2;
    distance_2 = lon2 - lon1;
    distance_1 = MAX(distance_1,distance_2);

    distance_2 = 360. - distance_1;

    results = MIN(distance_1,distance_2);

    return results;
}

int
check_if_in_box(float32 lat, float32 lon, float32 northern_boundary, float32 southern_boundary,
                float32 western_boundary, float32 eastern_boundary)
{
    int results;
    int east_west_results;
    int north_south_results;

    east_west_results = check_if_in_west_east_boundaries(lon, western_boundary, eastern_boundary);

    if (lat >= southern_boundary && lat <= northern_boundary)
    {
        north_south_results = TRUE;
    }
    else
    {
        north_south_results = FALSE;
    }

    if (north_south_results == TRUE && east_west_results == TRUE)
    {
        results = TRUE;
    }
    else 
    {
        results = FALSE;
    }

    return results;
}





int
check_if_in_west_east_boundaries(float32 lon, float32 western_boundary, float32 eastern_boundary)
{

    int results = FALSE;

    if (eastern_boundary >= western_boundary)
    {
        /* no date line crossing */

        if (lon >= western_boundary && 
            lon <= eastern_boundary) 
        {
            results = TRUE;
        }
    }
    else
    {
        /* date line crossing */

        if (lon >= western_boundary ||
            lon <= eastern_boundary)
        {
            results = TRUE;
        }
    }

    return results;
}






void
get_box_square_km(double northern_boundary, 
                     double southern_boundary, 
                     double western_boundary, 
                     double eastern_boundary,
                     double * box_square_km)
{

    double northern_box_square_km;
    double southern_box_square_km;


    if (northern_boundary > 0. && southern_boundary < .0)
    {
        /*************************************************************************************
        *    Divide into two boxes since equator is crossed
        *************************************************************************************/

        get_box_square_km_sublevel(northern_boundary,
                             (double) 0.,
                             western_boundary,
                             eastern_boundary,
                             &northern_box_square_km);

        get_box_square_km_sublevel((double) 0.,
                             southern_boundary,
                             western_boundary,
                             eastern_boundary,
                             &southern_box_square_km);

        *box_square_km = northern_box_square_km + southern_box_square_km;
    }
    else
    {
        get_box_square_km_sublevel(northern_boundary,
                             southern_boundary,
                             western_boundary,
                             eastern_boundary,
                             box_square_km);
    }
    
}

void
get_box_square_km_sublevel(double northern_boundary, 
                              double southern_boundary, 
                              double western_boundary, 
                              double eastern_boundary,
                              double * box_square_km)
{
    double lat_width_degrees;
    double lon_width_degrees;
    double north_lon_radius;
    double south_lon_radius;
    double north_lon_width_km;
    double south_lon_width_km;
    double avg_lon_width_km;
    double lat_width_km;

    if (eastern_boundary >= western_boundary)
    {
        /* no date line crossing */
        lon_width_degrees = eastern_boundary - western_boundary;
    }
    else
    {
        /* date line crossing */
        lon_width_degrees = 360. + eastern_boundary - western_boundary;
    }
 
    north_lon_radius = EARTH_RADIUS_EQUATORIAL * fabs( cos(northern_boundary * (PI/180.)) ); 
    north_lon_width_km = (PI/180.) * lon_width_degrees * north_lon_radius;

    south_lon_radius = EARTH_RADIUS_EQUATORIAL * fabs( cos(southern_boundary * (PI/180.)) ); 
    south_lon_width_km = (PI/180.) * lon_width_degrees * south_lon_radius;

    avg_lon_width_km = (north_lon_width_km + south_lon_width_km) / 2.;

    lat_width_degrees = northern_boundary - southern_boundary;
    lat_width_km = (PI/180.) * lat_width_degrees * EARTH_RADIUS_EQUATORIAL;

    *box_square_km = lat_width_km * avg_lon_width_km;

}




void
set_west_east_boundaries(float32 lon, float32 * western_boundary, float32 * eastern_boundary)
{
    float32 boundary_width;
    float32 width_to_west;
    float32 width_to_east;

    if (*western_boundary != -180. || *eastern_boundary != 180.) 
    {
        if (*eastern_boundary != DEFAULT_COORD_VALUE && *western_boundary != DEFAULT_COORD_VALUE)
        {
            if (check_if_in_west_east_boundaries(lon,*western_boundary,*eastern_boundary) != TRUE)
            {
                /************************************************************************
                *    Determine longitude width to east of current boundary
                ************************************************************************/
     
                if (lon > *eastern_boundary)
                {
                    /* dateline not crossed */
                    width_to_east = lon - *eastern_boundary;
                }
                else
                {
                    /* dateline crossed */
                    width_to_east = 360. + lon - *eastern_boundary;
                }
    
                /************************************************************************
                *    Determine longitude width to west of current boundary
                ************************************************************************/
     
                if (*western_boundary > lon)
                {
                    /* dateline not crossed */
                    width_to_west = *western_boundary - lon;
                }
                else
                {
                    /* dateline crossed */
                    width_to_west = *western_boundary + 360. - lon;
                }
    
                /************************************************************************
                *    Set closest west-east boundary
                ************************************************************************/
     
                if (fabs(width_to_west) <= fabs(width_to_east))
                {
                    *western_boundary = lon;
                }
                else
                {
                    *eastern_boundary = lon;
                }

                /************************************************************************
                *    Determine longitude width between western_boundary and eastern_boundary
                ************************************************************************/
     
                if (*eastern_boundary >= *western_boundary)
                {
                    /* no date line crossing */
                    boundary_width = *eastern_boundary - *western_boundary;
                }
                else
                {
                    /* date line crossing */
                    boundary_width = 360. + *eastern_boundary - *western_boundary;
                }

                /************************************************************************
                *    if west-to-east span > 355 then just set span to -180 to 180
                ************************************************************************/
     
                if (boundary_width > 355.)
                {
                    *western_boundary = -180.;
                    *eastern_boundary = 180.;
                }
            }
        }
        else
        {
            *eastern_boundary = lon;
            *western_boundary = lon;
        }
    }
}





